Efficient conjugate gradient based channel estimator

ABSTRACT

In forming a current channel estimate from a received signal y, the received signal y is decoded to form data s, a convolution matrix Ŝis formed from a first portion of the data s, a matrix F 1  is formed from a second portion the data s such that the second portion of the data s includes data that is less recent that the data in the first portion of the data s, a matrix F 2  is formed from a third portion the data s such that the third portion of the data s includes data that is more recent than the data in the first portion of the data s, a predicted channel estimate h pred  is determined, and a conjugate gradient algorithm is performed to determine the current channel estimate. The conjugate gradient algorithm is based on the received signal y, the matrix Ŝ, the matrices F 1  and F 2 , the predicted channel estimate h pred , and a previous channel estimate h 1 .

RELATED APPLICATIONS

This application is a continuation-in-part of U.S. application Ser. No.10/729,722 filed on Dec. 5, 2003 now U.S. Pat. No. 7,327,810.

TECHNICAL FIELD OF THE INVENTION

The present invention relates to a conjugate gradient based channelestimator.

BACKGROUND OF THE INVENTION

Since the adoption of the ATSC digital television (DTV) standard in 1996in the United States, there has been an ongoing effort to improve thedesign of receivers built for the ATSC DTV signal. The primary obstaclethat faces designers in designing receivers so that they achieve goodreception is the presence of multipath interference in the channel.

The broadcast television channel is a relatively severe multipathenvironment due to a variety of conditions that are encountered in thechannel and at the receiver. Channels are characterized by impulseresponses which may be several hundreds of symbols long, so that stronginterfering signals may arrive at the receiver both before and after thelargest amplitude signal. In addition, the signal transmitted throughthe channel is subject to time varying channel conditions due to themovement of the transmitter and signal reflectors, airplane flutter,and, for indoor reception, people walking around the room. If mobilereception is desired, movement of the receiver must also be considered.

Moreover, the ATSC DTV signal uses trellis coded 8-level vestigialsideband (usually referred to as 8T-VSB or, more simply, as 8-VSB) asthe modulation method. 8-VSB data symbols are real and have a signalpulse shape that is complex. Only the real part of the complex pulseshape is a Nyquist pulse. Therefore, even if there is no multipath, theimaginary part of the complex pulse shape contributes intersymbolinterference (ISI) when the channel gain seen by the equalizer is notreal.

Multipath and intersymbol interference adversely affect the ability ofthe receiver to correctly receive the symbols transmitted by thetransmitter. Therefore, designers add equalizers to receivers in orderto cancel the effects of multipath and intersymbol interference andthereby improve signal reception.

Because the channel is not known a priori at the receiver, the equalizermust be able to modify its response to match the channel conditions thatit encounters and to adapt to changes in those channel conditions. Toaid in the convergence of an adaptive equalizer to channel conditions,the field sync segment of the frame as defined in the ATSC standard maybe used as a training sequence for the equalizer. But when equalizationis done in the time domain, long equalizers (those having many taps) arerequired due to the long channel impulse responses that characterize thechannel.

The original Grand Alliance receiver used an adaptive decision feedbackequalizer (DFE) with 256 taps. This adaptive decision feedback equalizerwas adapted to the channel using a standard least mean square (LMS)algorithm, and was trained with the field sync segment of thetransmitted frame. The LMS algorithm converges quite slowly and, evenwith only 256 taps, does not always converge during a single trainingsequence. Because the field sync segment is transmitted relativelyinfrequently (about every 260,000 symbols), the total convergence timeof this equalizer is quite long if the equalizer only adapts on trainingsymbols prior to convergence.

Therefore, in order to adapt equalizers to follow channel variationsthat occur between training sequences, the addition of blind anddecision directed methods to equalizers has been suggested. However,when implemented in a realistic system, these methods may requireseveral data fields to achieve convergence, and convergence may not beachieved at all under difficult multipath conditions.

In any event, because multipath signals in the broadcast channel mayarrive many symbols after the main signal, the decision feedbackequalizer is invariably used in 8-VSB applications.

It has been argued that a blind decision feedback equalizer is requireddue to the rise in mean square error (MSE) between training sequences.However, adaptation of the trained equalizer in the simulation thatsupported this argument was frozen between training sequences. It ispossible that a decision-directed equalizer with good trackingperformance may be able to follow the channel variations tested.

Blind algorithms based on the Sato algorithm and on Godard's constantmodulus algorithm (CMA) have been proposed. The error term in both ofthese algorithms uses a continuous blending of a decision-directed termwith the blind term. This blending enables a smooth transition betweenthe blind mode and the decision-directed mode. However, when implementedin a realistic system, these algorithms also may take several datafields to converge.

As mentioned previously, adaptive equalizers utilizing the least meansquare (LMS) algorithm may converge slowly or not at all depending onthe channel conditions. Convergence may be adversely affected if theinput data auto-correlation matrix has a large eigenvalue spread. Also,if the decision feedback equalizer has not converged before the end ofthe training sequence, the shape of the objective function may change sothat it includes local minima. These local minima may be caused byclosed eye channel conditions and decision feedback equalizer errorpropagation.

The recursive least square (RLS) algorithm is well known to avoid theseconvergence problems. However, the recursive least square algorithm, inits basic form, requires the computationally intensive inversion of theinput data auto-correlation matrix.

Lattice based forms of the recursive least square algorithm avoid theneed for this matrix inversion. However, the lattice based forms of therecursive least square algorithm are not easily amenable to theadvantage of initialization from an initial channel impulse response(CIR) estimate, even though such an initialization may be desirable.

Recent work in reduced rank filtering has connected the multi-stagenested Wiener filter (MSNWF) of Goldstein and Reed to the conjugategradient algorithm (CG). It has been shown that the multi-stage nestedWiener filter solves the Wiener-Hopf equations in the Krylov subspaceassociated with the auto-correlation matrix and the cross-correlationvector. The multi-stage nested Wiener filter is then re-formulated usingthe Lanczos iteration. It has also been shown that the Lanczos-basedmulti-stage nested Wiener filter is equivalent to the conjugate gradientalgorithm. Since the multi-stage nested Wiener filter often needs fewdimensions to approach the performance of the full-rank Wiener filter,the conjugate gradient algorithm is a good candidate for an adaptiveequalization algorithm with fast convergence.

A description of the conjugate gradient optimization algorithm and manyof its mathematical properties may be found in “An Introduction toOptimization,” by E. K. P. Chong and S. Zak, New York, N.Y., John Wiley& Sons, 1996. The algorithm described is applicable to the optimizationof a fixed objective function. An excellent review of adaptive filteringwith the conjugate gradient algorithm is found in “Analysis of conjugategradient algorithms for adaptive filtering,” by P. S. Chang and A. N.Willson, Jr., vol. 48, pp. 409-418, IEEE Transactions on SignalProcessing, February 2000. The mathematical richness of the algorithmrelationships leads to a variety of options in the implementation of thebasic conjugate gradient algorithm for minimizing a quadratic objectivefunction. These options can be carried over to the adaptive situation toprovide a number of options for implementation of the algorithm. Onecharacteristic of the algorithm which may also be exploited is that itworks directly with the correlation matrices and vectors. Thischaracteristic makes initialization based on an initial channel impulseresponse straightforward, assuming information is available for thispurpose.

Like lattice based forms of the recursive least square algorithm, theconjugate gradient algorithm does not require the inversion of the inputdata auto-correlation matrix. Another characteristic of the conjugategradient algorithm may be exploited by recognizing that the input dataauto-correlation matrix is the product of two Toeplitz matrices.

It has been known internally by assignee to use these characteristics tosubstantially reduce the computation load in the conjugate gradientalgorithm and to even eliminate the need to form the inputauto-correlation matrix. Accordingly, the improved conjugate gradientalgorithm disclosed in that application may be used to achievesatisfactory convergence times for equalizers.

It has also been known internally by assignee that the conjugategradient algorithm can be used to estimate channels. The presentinvention reduces the computation complexity in performing the conjugategradient algorithm with respect to channel estimators.

SUMMARY OF THE INVENTION

In accordance with one aspect of the present invention, a method ofprocessing a received signal y to produce a current channel estimatefrom one or more past channel estimates that include a previous channelestimate h₁ comprises the following: (a) decoding the received signal yto form data s; (b) forming a convolution matrix Ŝ from a first portionof the data s; (c) forming a matrix F₁ from a second portion the data s,wherein the second portion of the data s includes data that is lessrecent than the data in the first portion of the data s; (d) forming amatrix F₂ from a third portion of the data s, wherein the third portionof the data s includes data that is more recent than the data in thefirst portion of the data s; (e) determining a predicted channelestimate h_(pred) based on the one or more past channel estimates; and,(f) performing a conjugate gradient algorithm to determine the currentchannel estimate, wherein the conjugate gradient algorithm is based onthe received signal y, the matrix Ŝ, the matrices F₁ and F₂, thepredicted channel estimate h_(pred), and the previous channel estimateh₁.

In accordance with another aspect of the present invention, a method ofprocessing a received signal y to produce a current channel estimatefrom one or more past channel estimates that include a previous channelestimate h₁ comprises the following: (a) decoding the received signal yto form data s; (b) forming a convolution matrix Ŝ from a first portionof the data s; (c) forming a matrix F₁ from a second portion of the datas, wherein the second portion of the data s includes data that is lessrecent than the data in the first portion of the data s; (d) forming amatrix F₂ from a third portion of the data s, wherein the third portionof the data s includes data that is more recent than the data in thefirst portion of the data s; (e) determining a predicted channelestimate h_(pred) based on the one or more past channel estimates; and,(f) performing a conjugate gradient algorithm to determine the currentchannel estimate, wherein the conjugate gradient algorithm is based onthe received signal y, the matrix Ŝ, the matrices F₁ and F₂, thepredicted channel estimate h_(pred), and the previous channel estimateh₁, wherein the conjugate gradient algorithm includes (i) forming FFTsbased on the received signal y, the matrix Ŝ, and the matrices F1 andF₂, (ii) multiplying the FFTs to form a multiplication product, and(iii) forming an inverse FFT of the multiplication product.

BRIEF DESCRIPTION OF THE DRAWINGS

These and other features and advantages will become more apparent from adetailed consideration of the invention when taken in conjunction withthe drawings in which:

FIG. 1 illustrates a filter that can be used in a decision feedbackequalizer according to an embodiment of the present invention;

FIG. 2 illustrates the decision feedback equalizer according to anembodiment of the present invention;

FIG. 3 illustrates a matrix S which is derived from a vector of symboldecisions of length 2496 and is useful in describing the presentinvention;

FIG. 4 illustrates the decomposition of the matrix S to form twomatrices Ŝ and F useful in implementing the present invention;

FIG. 5 shows a system for determining a channel estimate h based on thematrices Ŝ and F; and,

FIG. 6 illustrates the decomposition of the matrix F to form twomatrices F₁ and F₂.

DETAILED DESCRIPTION

As shown in FIG. 1, a filter 10 can be used as an equalizer and includesstorage registers 12[n], 12[n+N−3], 12[n+N−2], 12[n+N−1]. For example, Nmay be 512. The storage registers 12[n], . . . , 12[n+N−3], 12[n+N−2],12[n+N−1] store corresponding received data y[n], y[n+N−3], y[n+N−2],y[n+N−1]. The output of each of the storage registers 12[n], . . . ,12[n+N−3], 12[n+N−2], 12[n+N−1] is coupled to a corresponding multiplier14[n], 14[n+N−3], 14[n+N−2], 14[n+N−1]. Accordingly, the multipliers14[n], . . . , 14[n+N−3], 14[n+N−2], 14[n+N−1] receive the stored datafrom the storage registers 12[n], 12[n+N−3], 12[n+N−2], 12[n+N−1]. Themultipliers 14[n], . . . , 14[n+N−3], 14[n+N−2], 14[n+N−1] also receivecorresponding tap weights g[n], . . . , g[n+N−3], g[n+N−2], g[n+N−1].The multipliers 14[n], . . . , 14[n+N−3], 14[n+N−2], 14[n+N−1] multiplyeach of the stored data y[n], . . . , y[n+N−3], y[n+N−2], y[n+N−1] by acorresponding one of the tap weights g[n], . . . , g[n+N−3], g[n+N−2],g[n+N−1], and the multiplication results are summed by a summer 16 toprovide the output of the filter 10. Various algorithms as describedabove have been implemented in order to set the values of the tapweights g[n], . . . , g[n+N−3], g[n+N−2], g[n+N−1].

A decision feedback equalizer 20 is shown in FIG. 2 and includes a feedforward filter 22, a feedback filter 24, a summer 26, a hard decisiondevice 28, and a controller 30 that computes the tap weights for thefeed forward filter 22 and the feedback filter 24. Each of the feedforward filter 22 and the feedback filter 24 may comprise the filter 10of FIG. 1. The hard decision device 28, for example, may be a decoder ora slicer.

The output ŝ[n] of the summer 26 is the equalized output of the decisionfeedback equalizer 20. The hard decision device 28 (such as a slicer)determines the closest value {tilde over (s)}[n] for each of the symbolsprovided as the output ŝ[n] of the summer 26. These values {tilde over(s)}[n] are fed back as inputs to the feedback filter 24. The feedbackfilter 24 applies tap weights ĝ_(B)[n] to the values {tilde over (s)}[n]and provides its output to the summer 26. Likewise, the feed forwardfilter 22 applies tap weights ĝ_(F)[n] to the values y[n] of thereceived signal and provides its output to the summer 26. The summer 26sums the output of the feedback filter 24 with the output of the feedforward filter 22 in order to produce the output ŝ[n].

The output ŝ[n] of the decision feedback equalizer 20 is given by thefollowing equation:

$\begin{matrix}\begin{matrix}{{\hat{s}\lbrack n\rbrack} = {{{Re}\left\{ {g_{F}^{H}{y\lbrack n\rbrack}} \right\}} + {g_{B}^{T}{\overset{\sim}{s}\lbrack n\rbrack}}}} \\{= {{g_{FR}^{T}{y_{R}\lbrack n\rbrack}} + {g_{FI}^{T}{y_{I}\lbrack n\rbrack}} + {g_{B}^{T}{\overset{\sim}{s}\lbrack n\rbrack}}}}\end{matrix} & (1)\end{matrix}$where {tilde over (s)}[n] as indicated above is defined to be the symbolestimates provided by the hard decision device 28, whereg_(FR)=Re{g_(F)} is the real part of the tap weights g_(F) of the feedforward filter 22, where g_(FI)=Im{g_(F)} is the imaginary part of thetap weights g_(F) of the feed forward filter 22, and where g_(B) are thetap weights of the feedback filter 24.

The definitions given by the following equations may be assumed:{tilde over (y)} _(F) [n]=└{tilde over (y)} _(R) ^(T) [n]y _(I) ^(T)[n]┘ ^(T)  (2){tilde over (y)}[n]=└{tilde over (y)} _(F) ^(T) [n]{tilde over (s)} ^(T)[n]┘ ^(T)  (3){tilde over (g)} _(F) [n]=└g _(FR) ^(T) [n]g _(FI) ^(T) [n]┘ ^(T)  (4){tilde over (g)}[n]=└{tilde over (g)} _(F) ^(T) [n]{tilde over (g)} _(B)^(T) [n]┘ ^(T)  (5)where y_(R)=Re{y_(F)} is the real part of the received symbols y_(F),and where y₁=Im{y_(F)} is the imaginary part of the received symbolsy_(F). Then, the equalized symbols ŝ[n] provided by the decisionfeedback equalizer 20 are given by the following equation:

$\begin{matrix}\begin{matrix}{{\hat{s}\lbrack n\rbrack} = {{{\overset{\sim}{g}}_{F}^{T}{{\overset{\sim}{y}}_{F}\lbrack n\rbrack}} + {g_{B}^{T}{\overset{\sim}{s}\lbrack n\rbrack}}}} \\{= {{\overset{\sim}{g}}^{T}{\overset{\sim}{y}\lbrack n\rbrack}}}\end{matrix} & (6)\end{matrix}$The equalizer input y[n] and the feed forward tap weights g_(F) arecomplex values, the feedback tap weights g_(B) are real, and the outputŝ[n] of the decision feedback equalizer 20 is real. The imaginary partof the output of the feed forward filter 22, g_(FI), may be ignored and,therefore, not computed.

It is desired to determine the tap weight vector {tilde over (g)}[n] sothat the error at the output of the decision feedback equalizer 20 isminimized according to the following equations:

$\begin{matrix}{{\overset{\sim}{g}\lbrack n\rbrack} = {\begin{matrix}{\arg\mspace{11mu}\min} \\g\end{matrix}\frac{1}{2}E\left\{ \left( {{s\lbrack n\rbrack} - {\hat{s}\lbrack n\rbrack}} \right)^{2} \right\}}} & (7) \\\begin{matrix}{{F(g)} = {E\left\{ \left( {{s\lbrack n\rbrack} - {\hat{s}\lbrack n\rbrack}} \right)^{2} \right\}}} \\{= {\sigma_{s}^{2} - {{{\overset{\sim}{g}}^{H}\lbrack n\rbrack}r_{s\overset{\sim}{y}}} - {r_{s\overset{\sim}{y}}^{H}{\overset{\sim}{g}\lbrack n\rbrack}} + {{{\overset{\sim}{g}}^{H}\lbrack n\rbrack}R_{\overset{\sim}{y}\overset{\sim}{y}}{g\lbrack n\rbrack}}}}\end{matrix} & (8)\end{matrix}$where σ_(s) ² is the symbol variance and is equal to the average symbolenergy, R_({tilde over (y)}{tilde over (y)}) is an autocorrelationmatrix, and r_(s{tilde over (y)}) is a cross correlation vector. Theautocorrelation matrix R_({tilde over (y)}{tilde over (y)}) and thecross correlation vector r_(s{tilde over (y)}) are defined in detailhereinafter. The error function F(g) is minimized by setting itsgradient to zero according to the following equation:∇F(g)=R _({tilde over (y)}{tilde over (y)}) {tilde over (g)}[n]− r_(s{tilde over (y)})=0  (9)This minimization leads to the following well known Wiener-Hopfequation:R _({tilde over (y)}{tilde over (y)}) {tilde over (g)}[n]= r_(s{tilde over (y)})  (10)Solving for the desired tap weights {tilde over (g)}[n] leads to thefollowing equation:

$\begin{matrix}{{\overset{\sim}{g}\lbrack n\rbrack} = {R_{\overset{\sim}{y}\overset{\sim}{y}}^{- 1}r_{s\overset{\sim}{y}}}} & (11)\end{matrix}$So, in order to calculate the desired tap weight vector, the largesystem of equations represented by equation (11) must be solved. Thissolution requires the inversion of the large matrixR_({tilde over (y)}{tilde over (y)}) and, therefore, a considerableamount of computation.

One method for reducing the amount of computation involved in solving alarge system of equations is the aforementioned conjugate gradientalgorithm. In order to set up equation (11), estimation of theauto-correlation matrix R_(yy)[n=0] and the cross correlation vectorr_(sy)[n=0] is required. By definition, the auto-correlation matrixR_(yy)[n] is given by the following equation:R _({tilde over (y)}{tilde over (y)}) [n]E[{tilde over (y)}[n]{tildeover (y)} ^(H) [n]]  (12)and the cross correlation vector r_(sy)[n] is given by the followingequation:r _(s{tilde over (y)}) [n]=E[s[n]{tilde over (y)}[n]]  (13)where s[n] represents the transmitted symbols, E is the expectationoperator, and H denotes the Hermitian transpose.

The auto-correlation matrix R_({tilde over (y)}{tilde over (y)})[n] andthe cross correlation vector r_(s{tilde over (y)})[n] can be estimatedusing the forgetting factor method according to the following equations:

$\begin{matrix}\begin{matrix}{{{\hat{R}}_{\overset{\sim}{y}\overset{\sim}{y}}\lbrack n\rbrack} = {\sum\limits_{m = 0}^{n - 1}{\gamma^{m}{\overset{\sim}{y}\left\lbrack {n - m} \right\rbrack}{{\overset{\sim}{y}}^{H}\left\lbrack {n - m} \right\rbrack}}}} \\{= {{\gamma{{\hat{R}}_{yy}\left\lbrack {n - 1} \right\rbrack}} + {{\overset{\sim}{y}\lbrack n\rbrack}{{\overset{\sim}{y}}^{H}\lbrack n\rbrack}}}}\end{matrix} & (14) \\\begin{matrix}{{{\hat{r}}_{s\overset{\sim}{y}}\lbrack n\rbrack} = {\sum\limits_{m = 0}^{n - 1}{\gamma^{m}{s\left\lbrack {n - m} \right\rbrack}{{\overset{\sim}{y}}^{H}\left\lbrack {n - m} \right\rbrack}}}} \\{= {{\gamma{{\hat{r}}_{s\overset{\sim}{y}}\left\lbrack {n - 1} \right\rbrack}} + {{s\lbrack n\rbrack}{{\overset{\sim}{y}}^{H}\lbrack n\rbrack}}}}\end{matrix} & (15)\end{matrix}$where 0<γ<1 is the forgetting factor. If γ is set to 0, the equations(14) and (15) yield the estimates of the auto-correlation matrix{circumflex over (R)}_({tilde over (y)}{tilde over (y)}) and the crosscorrelation vector {circumflex over (r)}_(s{tilde over (y)}) from themost recent sample vector y.

In determining the tap weights g_(F) for the feed forward filter 22 andg_(B) for the feedback filter 24, the controller 30 initializes theauto-correlation matrix R_({tilde over (y)}{tilde over (y)})[0] usingequation (14) and the cross correlation vector r_(s{tilde over (y)})[0]using equation (15) by performing averaging over a period of time suchas ½ data field or so. The controller 30 also initializes the tapweights {tilde over (g)}₀[0]=0 and initializes the decision vector s[0]to the first N_(fb) of the known training symbols that are contained inthe field sync portion of each data field as discussed above.

Further, the controller 30 initializes the data vector {tilde over(y)}[0] as follows. The data vector {tilde over (y)}_(F)[0] as given inequation (2) is initialized by setting the real party_(R)[n=0]={y_(R)[0], . . . , y_(R)[N_(ff)−1]}^(T) such that y_(R)[0], .. . , y_(R)[N_(ts)−N_(fb)−1] are the real part of the remaining trainingdata (after those in {tilde over (s)}[0]) and y_(R)[N_(ts)−N_(fb)], . .. , y_(R)[N_(ff)−1] are the real part of the received data after thetraining sequence. This setting is repeated for {tilde over (y)}_(I)[0].Then, {tilde over (y)}_(F)[0]=[y_(R) ^(T)[0]y_(I) ^(T)[0]]^(T) isconstructed as given in equation (2) and {tilde over (y)}[0]=[{tildeover (y)}_(F) ^(T)[0]{tilde over (s)}^(T)[0]]^(T) is constructed asgiven in equation (3) based on the initialized vectors {tilde over(y)}_(F)[0], {tilde over (y)}_(R)[0], and {tilde over (y)}_(I)[0] asgiven above. In an ATSC digital television environment, for example, thenumber of training symbols N_(ts) may be 728, the number of feed forwardtaps N_(ff) may be 512, and the number of feedback taps N_(fb) may alsobe 512.

Further, the controller 30 initializes the output ŝ[n] of the decisionfeedback equalizer 20 as ŝ[0]={tilde over (g)}₀ ^(T)[0]{tilde over(y)}[0] in accordance with equation (6), the controller 30 initializes agradient t_(n)[n] as t₀[0]={circumflex over(R)}_({tilde over (y)}{tilde over (y)})[0]{tilde over(g)}₀[0]−{circumflex over (r)} _(sy)[0], and the controller 30 sets thevariable n to 0.

The controller 30 then executes a conjugate gradient algorithm once foreach received data element, where each iteration of the algorithm isinitiated upon the shifting of a new received data element into the feedforward filter 22. The parameters of the conjugate gradient algorithm asdescribed above are defined in the table.

TABLE Parameter Name Equalizer Notation Object Function Argumentg_(k)[n] (Equalizer) Argument Step Size α_(k)[n] Gradient t_(k)[n]Conjugate Directions b_(k)[n] Conjugate Direction Step β_(k)[n] Size

This conjugate gradient algorithm is as follows:

-   -   1. calculate conjugate direction b_(k)[n]        -   if n modulo(N_(taps))=0            -   b₀[n]=−t₀[0] (reset)        -   else

${\beta_{0}\lbrack n\rbrack} = \frac{\left( {{t_{o}\lbrack n\rbrack} - {t_{- 1}\lbrack n\rbrack}} \right)^{H}{t_{0}\lbrack n\rbrack}}{\left( {{t_{0}\lbrack n\rbrack} - {t_{- 1}\lbrack n\rbrack}} \right)^{H}{b_{- 1}\lbrack n\rbrack}}$

-   -   -   b₀[n]=−t₀[n]+β₀[n]b⁻¹[n]

    -   2. calculate argument step size α_(k)[n]

${\alpha_{0}\lbrack n\rbrack} = \frac{{- {b_{0}^{H}\lbrack n\rbrack}}{t_{0}\lbrack n\rbrack}}{{b_{0}^{H}\lbrack n\rbrack}{{\hat{R}}_{\overset{\sim}{y}\overset{\sim}{y}}\lbrack n\rbrack}{b_{0}\lbrack n\rbrack}}$

-   -   3. tap update        -   {tilde over (g)}[n]={tilde over (g)}₀[n]+α₀[n]b₀[n]    -   4. n=n+1        -   {tilde over (g)}₀[n]={tilde over (g)}₁[n−1]        -   t⁻¹[n]=t₀[n−1]        -   b⁻¹[n]=b₀[n−1]    -   5. update input vector {tilde over (y)}[n] with the next        received data element and update {circumflex over        (R)}_({tilde over (y)}{tilde over (y)})[n] per equation (14)    -   6. Calculate output of the decision feedback equalizer 20        -   ŝ[n]={tilde over (g)}₀ ^(T)[n]{tilde over (y)}[n]    -   7. calculate gradient t_(k)[n]        -   t₀[n]=t⁻¹[n]+α₀[n]└{circumflex over            (R)}_({tilde over (y)}{tilde over (y)})[n]b₀[n]    -   8. Return to step 1

The largest of the part of the computations required by this conjugategradient algorithm is contributed by the matrix-vector multiplications{circumflex over (R)}_({tilde over (y)}{tilde over (y)})[n]b₀[n] ofsteps 2 and 7. These multiplications require the update of matrix{circumflex over (R)}_({tilde over (y)}{tilde over (y)})[n] from thevector {tilde over (y)}[n] as in equation (14). This update requires onemultiply for each of the N_(taps) ² matrix elements, a total of N_(taps)² multiplies. Then, the matrix-vector multiplication requires anadditional N_(taps) ² multiplies. Accordingly, the total computation isN_(taps) ²+N_(taps) ²=2N_(taps) ² multiplies.

This computational load can be decreased as explained below.

A matrix-vector multiplication given by a=Yb, where Y is a convolutionmatrix, b is a vector, and a is the multiplication result, is equivalentto a[n]=y[n]*b[n] where * denotes linear convolution. These equationsmay be represented in matrix form by the following equation:

$\begin{matrix}{\begin{bmatrix}{a\lbrack 0\rbrack} \\{a\lbrack 1\rbrack} \\{a\lbrack 2\rbrack}\end{bmatrix} = {\begin{bmatrix}{y\lbrack 0\rbrack} & {y\left\lbrack {- 1} \right\rbrack} & {y\left\lbrack {- 2} \right\rbrack} \\{y\lbrack 1\rbrack} & {y\lbrack 0\rbrack} & {y\left\lbrack {- 1} \right\rbrack} \\{y\lbrack 2\rbrack} & {y\lbrack 1\rbrack} & {y\lbrack 0\rbrack}\end{bmatrix}\begin{bmatrix}{b\lbrack 0\rbrack} \\{b\lbrack 1\rbrack} \\{b\lbrack 2\rbrack}\end{bmatrix}}} & (16)\end{matrix}$where n is chosen as a small number for explanatory purposes only. Thecolumns and rows of the matrix Y in equation (16) can be circularlyextended in accordance with the following equation:

$\begin{matrix}{\begin{bmatrix}{a\lbrack 0\rbrack} \\{a\lbrack 1\rbrack} \\{a\lbrack 2\rbrack} \\{dc} \\{dc}\end{bmatrix} = {\begin{bmatrix}{y\lbrack 0\rbrack} & {y\left\lbrack {- 1} \right\rbrack} & {y\left\lbrack {- 2} \right\rbrack} & {y\lbrack 2\rbrack} & {y\lbrack 1\rbrack} \\{y\lbrack 1\rbrack} & {y\lbrack 0\rbrack} & {y\left\lbrack {- 1} \right\rbrack} & {y\left\lbrack {- 2} \right\rbrack} & {y\lbrack 2\rbrack} \\{y\lbrack 2\rbrack} & {y\lbrack 1\rbrack} & {y\lbrack 0\rbrack} & {y\left\lbrack {- 1} \right\rbrack} & {y\left\lbrack {- 2} \right\rbrack} \\{y\left\lbrack {- 2} \right\rbrack} & {y\lbrack 2\rbrack} & {y\lbrack 1\rbrack} & {y\lbrack 0\rbrack} & {y\left\lbrack {- 1} \right\rbrack} \\{y\left\lbrack {- 1} \right\rbrack} & {y\left\lbrack {- 2} \right\rbrack} & {y\lbrack 2\rbrack} & {y\lbrack 1\rbrack} & {y\lbrack 0\rbrack}\end{bmatrix}\begin{bmatrix}{b\lbrack 0\rbrack} \\{b\lbrack 1\rbrack} \\{b\lbrack 2\rbrack} \\0 \\0\end{bmatrix}}} & (17)\end{matrix}$where dc means don't care. Equation (17) is equivalent to a[n]=y[n]{circle around (×)}b[n] where {circle around (×)} denotes circularconvolution. Circular convolution is equivalent to the followingoperation:γ=FFT(y)  (18)

=FFT(b)  (19)

=γ×

  (20)a=IFFT(

)  (21)where FFT in equations (18) and (19) indicates the Fast FourierTransform operation, IFFT in equation (21) indicates the Inverse FastFourier Transform operation, and x in equation (20) indicates point-wisevector multiplication such that the elements of the resultant vector aregiven by Y₀×B₀, Y₁×B₁, Y₂×B₂, etc.

An embellishment on the above process involves 0 padding the columns tolengths of 2^(n). This 0 padding, as is well known, permits acomputationally efficient FFT. Accordingly, the matrix-vectormultiplication a=Yb can be written with 0 padding in matrix formaccording to the following equation:

$\begin{matrix}{\begin{bmatrix}{a\lbrack 0\rbrack} \\{a\lbrack 1\rbrack} \\{a\lbrack 2\rbrack} \\{dc} \\{dc} \\{dc} \\{dc} \\{dc}\end{bmatrix} = {\begin{bmatrix}{y\lbrack 0\rbrack} & {y\left\lbrack {- 1} \right\rbrack} & {y\left\lbrack {- 2} \right\rbrack} & 0 & 0 & 0 & {y\lbrack 2\rbrack} & {y\lbrack 1\rbrack} \\{y\lbrack 1\rbrack} & {y\lbrack 0\rbrack} & {y\left\lbrack {- 1} \right\rbrack} & {y\left\lbrack {- 2} \right\rbrack} & 0 & 0 & 0 & {y\lbrack 2\rbrack} \\{y\lbrack 2\rbrack} & {y\lbrack 1\rbrack} & {y\lbrack 0\rbrack} & {y\left\lbrack {- 1} \right\rbrack} & {y\left\lbrack {- 2} \right\rbrack} & 0 & 0 & 0 \\0 & {y\lbrack 2\rbrack} & {y\lbrack 1\rbrack} & {y\lbrack 0\rbrack} & {y\left\lbrack {- 1} \right\rbrack} & {y\left\lbrack {- 2} \right\rbrack} & 0 & 0 \\0 & 0 & {y\lbrack 2\rbrack} & {y\lbrack 1\rbrack} & {y\lbrack 0\rbrack} & {y\left\lbrack {- 1} \right\rbrack} & {y\left\lbrack {- 2} \right\rbrack} & 0 \\0 & 0 & 0 & {y\lbrack 2\rbrack} & {y\lbrack 1\rbrack} & {y\lbrack 0\rbrack} & {y\left\lbrack {- 1} \right\rbrack} & {y\left\lbrack {- 2} \right\rbrack} \\{y\left\lbrack {- 2} \right\rbrack} & 0 & 0 & 0 & {y\lbrack 2\rbrack} & {y\lbrack 1\rbrack} & {y\lbrack 0\rbrack} & {y\left\lbrack {- 1} \right\rbrack} \\{y\left\lbrack {- 1} \right\rbrack} & {y\left\lbrack {- 2} \right\rbrack} & 0 & 0 & 0 & {y\lbrack 2\rbrack} & {y\lbrack 1\rbrack} & {y\lbrack 0\rbrack}\end{bmatrix}\begin{bmatrix}{b\lbrack 0\rbrack} \\{b\lbrack 1\rbrack} \\{b\lbrack 2\rbrack} \\0 \\0 \\0 \\0 \\0\end{bmatrix}}} & (22)\end{matrix}$It is noted that matrix Y has a Toeplitz structure, and that Y^(H) isalso a Toeplitz matrix. (H denotes Hermitian transpose.)

It will be shown that the need for both theR_({tilde over (y)}{tilde over (y)}) update and theR_({tilde over (y)}{tilde over (y)}) multiply can be eliminated, andthat the matrix R_({tilde over (y)}{tilde over (y)}) need not even becreated.

Let γ=0 for equation (14). It is well known thatR _({tilde over (y)}{tilde over (y)}) ={tilde over (y)}{tilde over (y)}^(H) =YY ^(H)where Y is a convolution matrix based on the vector y (that is, y formsthe first column of the matrix Y). The vector {tilde over (y)} may be 0padded to length 2^(n). Let {tilde over (y)}′ be the first column ofY^(H). For the conjugate gradient operation, the calculation given bythe following equation must be made:ν=R _({tilde over (y)}{tilde over (y)}) [n]b ₀ [n]=Y ^(H) Yb ₀ [n]  (23)In order to calculate ν, the quantity Y′b₀[n] should first be calculatedaccording to the following equations:γ=FFT(y)  (24)

=FFT(b)  (25)X=γ×

  (26)x=IDFT(X)  (27)x_(N)=first N elements of x, remaining elements=0  (28)X _(N) =FFT(x _(N))  (29)where x in equation (26) denotes point-wise vector multiplication. Then,Y^(H)×X_(N) is calculated according to the following equations:γ′=FFT(y′)  (30)

=γ′×X   (31)ν=IDFT

  (32)ν_(N)=first N elements of ν, remaining elements=0  (33)The number of multiplies in this process, particularly for largeN_(taps), is much less that 2N_(taps) ².

In the process given by equations (24)-(33), b is the vector b₀ of theeight step algorithm described above, y is a vector of received symbolsas described above, and ν_(N) is the quantity {circumflex over(R)}_({tilde over (y)}{tilde over (y)})[n]b₀[n] in steps 2 and 7 of theconjugate gradient algorithm described above.

The conjugate gradient algorithm described above can be used to makeother calculations. For example, the conjugate gradient algorithmdescribed above can be used to calculate updated channel impulseresponse estimates. This calculation can be useful, for instance,because known alternative methods of calculating MMSE (Minimum MeanSquare Error) equalizer tap weights are based on having a channelimpulse response estimate.

The conjugate gradient method described above permits the continuouscalculation of channel impulse response estimate updates betweentraining sequences (when only a priori unknown symbols are received) ina receiver for a time varying multipath channel. This calculationutilizes a least squares channel impulse response estimation calculatedfrom an over-determined system of equations. Dynamic changes to thechannel impulse response may be tracked by using receiver trellisdecoder decisions on input symbols to form a long sequence of almostperfectly known symbols. This sequence should have relatively fewerrors.

The channel impulse response to be estimated is of lengthL_(h)=L_(ha)+L_(hc)+1 where L_(ha) is the length of the anti-causal partof the channel impulse response and L_(hc) is the length of the causalpart of the channel impulse response.

The reliable trellis decoder decisions on the input symbols is a vectordesignated herein as s of length L_(s). A Toeplitz matrix S may bedefined based on this vector s according to the following equation:

$\begin{matrix}{S = \begin{bmatrix}S_{{Lh} - 1} & S_{{Lh} - 2} & — & — & — & — & — & S_{0} \\— & S_{{Lh} - 1} & — & — & — & — & — & — \\— & — & \; & \; & \; & \; & \; & — \\— & — & \; & \; & \; & \; & \; & — \\— & — & \; & \; & \; & \; & \; & S_{{Lh} - 1} \\— & — & \; & \; & \; & \; & \; & — \\— & — & \; & \; & \; & \; & S_{{Ls} - {Lh}} & — \\S_{{Ls} - 1} & S_{{Ls} - 2} & — & — & — & — & — & S_{{Ls} - {Lh}}\end{bmatrix}} & (34)\end{matrix}$where the elements in the matrix of equation (34) are real and consistof the symbol decisions of vector s.

To ensure an over-determined system of equations, the followinginequality is necessary:L _(S)>2L _(h)−1.  (35)The matrix S is of dimension (L_(S)−L_(h)+1)×L_(h) with(L_(S)−L_(h)+1)≧L_(h). The received signal vector is y with dataelements y_(i) for L_(hc)≦i≦(L_(S)−L_(ha)−1) where y_(i) is the receivedsymbol corresponding to input symbol decision s_(i). Then, y=Sh+u whereh is the L_(h) long channel impulse response vector and u is a whitenoise vector.

The least squares solution for h is given by the following equation:ĥ_(LS)=(S ^(T) S)⁻¹ S ^(T) y.  (36)

By utilizing reliable trellis decoder input symbol decisions, thereshould be sufficient support for calculating a channel impulse responseestimate of the required length. As required by inequality (35), thevector s of symbol decisions must be at least twice as long as thechannel impulse response being estimated.

For an 8 VSB receiver application, the following parameters may beassumed: L_(h)=512, L_(ha)=63, L_(hc)=448, and L_(S)=2496. The vector smay be formed from a sequence of trellis decoder decisions on the inputsymbols. Normally, the trellis decoder would just make output bit pairdecisions, but it can also make equally reliable decisions on the inputsymbols.

The vector s, for example, may be selected as 3 segments (L_(S)=2496symbols) long. So, three data segments may be used to produce a singlechannel impulse response estimate update. A new channel impulse responseestimate update can be obtained once per segment by proceeding in asliding window manner. Optionally, several consecutive channel impulseresponse estimate updates can be averaged in order to further improvechannel impulse response estimate accuracy if necessary. Of course, thisadditional step of averaging can be a problem if the channel impulseresponse is varying rapidly.

A vector b with fewer than 3 segments of symbol decisions may be used,but as stated in inequality (35), the length of vector b must be atleast twice as long as the channel impulse response to be estimated.Longer b vectors help to diminish the adverse effect of additive whiteGaussian noise on the channel.

The system of equations represented by equation (36) may be solved usingthe conjugate gradient algorithm in a manner similar to that previouslydescribed to solve equation (11) for the MMSE tap weights where S^(T)yin equation (36) takes the place of r_(s{tilde over (y)}) in equation(11), where S^(T)S in equation (36) takes the place ofR_({tilde over (y)}{tilde over (y)}) in equation (11), and where ĥ_(LS)in equation (36) takes the place of {tilde over (g)} in equation (11).

So, the conjugate gradient algorithm described above can facilitate thecomputation of (S^(T)S)b for the case of the channel impulse responseestimation just as it facilitated the computation ofR_({tilde over (y)}{tilde over (y)})b for the calculation of tapweights.

In further reducing the computation complexity necessary to implementthe matrix-vector multiplications required to perform the conjugategradient algorithm, such as when the conjugate gradient algorithm isused to estimate a channel as described above, it is useful to note thatthe conjugate gradient algorithm may be used to solve the followinglinear system for the channel estimate h:Sh=y.  (37)The matrix S is an m×n (e.g., 1985×512) Toeplitz matrix derived fromsymbol decisions, and y is a vector (the observation vector) derivedfrom the received data. The symbol decisions can be made, for example,by a trellis decoder or a slicer.

The system may be written as S^(T)Sh=S^(T)y by multiplying both sides ofequation (37) by S^(T). The conjugate gradient algorithm can be used tosolve this system for the channel estimate h according to the followingsteps:

-   (1) Calculate the residual vector r₁ representing the system    S^(T)y−S^(T)Sh₁ as follows:    r ₁ =S ^(T) y−S ^(T) Sh ₁  (38)    where r₁ is a 512×1 vector, and where h₁ is the previous channel    estimate which may be initialized to any desired value such as h₁=0.-   (2) For k=1 to n, iteratively calculate the following:    (a) d _(k) =r _(k)+β_(k) d _(k−1)  (39)    where β₁=0,

${\beta_{k \geq 2} = \frac{r_{k}^{T} \cdot r_{k}}{r_{k - 1}^{T} \cdot r_{k - 1}}},$and where d_(k) is a 512×1 vector;(b) h _(k+1) =h _(k) +α _(k) d _(k)  (40)where h_(k+1) is a 512×1 vector, where α_(k) is 1×1, where

${\alpha_{k} = \frac{r_{k}^{T} \cdot r_{k}}{d_{k} \cdot q_{k}}},$where q_(k)=S^(T)Sd_(k), where S^(T) is a 512×1985 matrix, and where Sis a 1985×512 matrix; and,(c) r _(k+1) =r _(k) −α _(k) q _(k)  (41)where r_(k+1) is a 512×1 vector.

Steps 2(a)-(c) are repeated enough times until a desired accuracy isachieved for the channel estimate h. The index k, for example, may beincremented up to a value of 5 although any other desired value can beused.

It will be seen that step 2(b) includes the matrix vector multiplicationS^(T)Sd_(k). This operation is computationally very complex. A techniquefor reducing the computational complexity of implementing this operationis suggested above with respect to equations (23)-(33).

A straightforward approach to this technique would be to form the FFT ofeach of the three components S^(T), S, and d_(k), point-wise vectormultiply the three resulting FFTs, and then take the inverse FFT of themultiplication result. This series of calculations, however, does notprovide an accurate answer because the matrix S is not a trueconvolution matrix, and the IFFT of the triple product, therefore,produces a corrupted result.

FIG. 3 illustrates the matrix S which is derived from a vector of symboldecisions of length 2496 (three segments of a VSB field). The first rowof the matrix S contains symbol decision 512 in the first column down tosymbol decision 1 in the last column, the second row of the matrix Scontains symbol decision 513 in the first column down to symbol decision2 in the last column, and so on. The portions of the matrix S whichprevent it from being a true convolution matrix are the triangular areas40 and 42. The diagonal lines in FIG. 3 are parallel.

According to the present invention, and as shown in FIG. 4, matrix S ofFIG. 3 is decomposed to form two matrices Ŝ and F such that the matrix Ŝis a true convolution matrix. Accordingly, the area 44 of the matrix Ŝis the same as the area 44 of the matrix S, and the symbol decisions inthe areas 40 and 42 of the matrix Ŝ are set to zero. The symboldecisions in the area 44 of the matrix F are set to zero, and the areas40 and 42 of the matrix F are the same as the areas 40 and 42 of thematrix S. A new linear system based on Ŝ and F is then defined which maybe solved for the channel estimate h. Advantageously, FFT operations ofŜ may be used for implementing the matrix vector multiplication of theconjugate gradient algorithm without corruption because Ŝ is a trueconvolution matrix.

More specifically, because the linear system is defined by y=Sh, thefollowing equations result:y=(Ŝ+F)h  (42)y=Ŝh+Fh  (43)y−Fh=Ŝh  (44)Defining y−Fh₁=ŷ, where h₁ is a previous channel estimate derived by anyknown technique, the new system expression becomes ŷ=Ŝh. This expressioncan be written as Ŝ^(T)ŷ=Ŝ^(T)Ŝh. The conjugate gradient algorithm canthen be applied to solve for h in the new system as follows:

-   -   (1) ŷ=y−Fh₁,        -   r₁=Ŝ^(T)ŷ−Ŝ^(T)Ŝh₁    -   (2) For k=1 to n, iteratively calculate        -   (a) d_(k)=r_(k)+β_(k)d_(k−1)        -   (b) h_(k+1)=h_(k)+α_(k)d_(k)        -   (c) r_(k+1)=r_(k)−α_(k)q_(k−1)            where h₁ may be initialized to any desired value such as 0,            and where d₀, r_(o), β_(k), and α_(k) are as previously            defined.

As indicated above following equation (40), step 2(b) includescalculating the expression q_(k)=Ŝ^(T)Ŝd_(k). This multiplication may beconveniently solved by forming the FFT of (i) the matrix Ŝ, (ii) d_(k),and (iii) Ŝ^(T). The three FFT results are then point-wise vectormultiplied, and the inverse FFT (IFFT) is taken of the multiplicationresult. The first N elements of this IFFT produces q_(k). N, forexample, can have a value of 512 for a 1024 point FFT/IFFT. However, Ncan have other values as desired.

FIG. 5 shows a system 50 for implementing this technique so as todetermine a channel estimate h. The system 50 includes a decoder 52 thatdecodes the received data y to produce the symbols s. As describedabove, the decoder 52 may be a trellis decoder. However, the decoder 52may be other forms of decoders such as slicers. A channel estimator 54produces a channel estimate h using the techniques described above.

As discussed above, the received data vector y may be related to thesymbol decisions s by the following equation:y=Sh+u  (45)where S is matrix formed of the symbol decisions and is of dimension(L_(S)−L_(h)+1)×L_(h) with (L_(S)−L_(h)+1)≧L_(h), where h is the L_(h)long channel impulse response vector, and where u is a white noisevector. Ignoring the noise vector, equation (45) may be rewritten asequations (42)-(44) where S is shown in FIG. 3, Ŝ and F are shown inFIG. 4, and S=Ŝ+F. As given above, ŷ=y−Fh₁ where h₁ is the previouschannel estimate, and the following new system of equations is solvedfor the channel estimate:ŷ=Ŝh  (46)

However, the previous channel estimate h₁ may be too outdated to be usedin the equation ŷ=y−Fh₁. More particularly, the lower left hand cornerof the matrix F has the newest (most recent) symbol decisions, and theprevious channel estimate h₁ may not be valid for these most recentsymbol decisions, such as in a dynamic channel environment wheremultipath taps may have Doppler shifts. It is, therefore, desirable tohave a more accurate channel vector than the previous estimate, whichmay be too old.

A simple remedy to this problem is to separate the matrix F into twomatrices as defined in the following equations:F=F ₁ +F ₂  (47)where FIG. 6 illustrates the decomposition of the matrix F to form thetwo matrices F₁ and F₂. In addition, a predicted channel vector denotedh_(pred) may be redefined according to the following equation:ŷ=y−F ₁ h ₁ −F ₂ h _(pred)  (48)Then, equation (46) can be solved for h using the new ŷ in the conjugategradient algorithm and FFT process described above.

Several methods can be used to determine the predicted channel estimateh_(pred). For example, a simple linear extrapolation may be used basedon the use of past k channel estimates where k≧2. If k=2, for example,the predicted channel vector can be obtained by defining a differencevector Δ according to the following equation:Δ=h _(k−1) −h _(k−2)  (49)Then, the following equations result:h _(pred) =h _(k−1)+Δ  (50)h _(pred) =h _(k−1)+(h _(k−1) −h _(k−2))  (51)h _(pred)=2h _(k−1) −h _(k−2)  (52)In the case where k=2, equation (52) requires the storage of k=2 pastchannel estimations, namely h_(k−1) and h_(k−2). It is noted thath₁=h_(k−1). Using simple linear extrapolation can be justified byconsidering a ghost (multipath) caused by a vehicle (truck, airplane,etc.) moving at a nearly constant speed. In this case, if the time spanbetween two consecutive channel estimates are close enough, a linearextrapolation will work satisfactorily. In the case of digitaltelevision, the channel estimate is updated once per 832 symbols (i.e.,once per segment or every 77 microseconds) so that a linearinterpolation should be sufficient to improve the channel estimate.

Modifications of the present invention will occur to those practicing inthe art of the present invention. For example, the Fast FourierTransform (FFT) and Inverse Fast Fourier Transform (IFFT) are used inthe present invention as described above. Alternatively, the DiscreteFourier Transform (DFT) and Inverse Discrete Fourier Transform (IDFT)can be used instead of the FFT and IFFT.

Accordingly, the description of the present invention is to be construedas illustrative only and is for the purpose of teaching those skilled inthe art the best mode of carrying out the invention. The details may bevaried substantially without departing from the spirit of the invention,and the exclusive use of all modifications which are within the scope ofthe appended claims is reserved.

1. A method of processing a received signal y to produce a currentchannel estimate from one or more past channel estimates that include aprevious channel estimate h₁ comprising: (a) decoding the receivedsignal y to form data s; (b) forming a convolution matrix Ŝ from a firstportion of the data s; (c) forming a matrix F₁ from a second portion thedata s, wherein the second portion of the data s includes data that isless recent than the data in the first portion of the data s; (d)forming a matrix F₂ from a third portion the data s, wherein the thirdportion of the data s includes data that is more recent than the data inthe first portion of the data s; (e) determining a predicted channelestimate h_(pred) based on the one or more past channel estimates; and,(f) performing a conjugate gradient algorithm to determine the currentchannel estimate, wherein the conjugate gradient algorithm is based onthe received signal y, the matrix Ŝ, the matrices F₁ and F₂, thepredicted channel estimate h_(pred), and the previous channel estimateh₁.
 2. The method of claim 1 wherein the determining of a predictedchannel estimate comprises extrapolating the predicted channel estimateh_(pred) from k of the past channel estimates.
 3. The method of claim 2wherein k≧2.
 4. The method of claim 1 wherein the received signal y, thematrix Ŝ, the matrices F₁ and F₂, the predicted channel estimateh_(pred), and the previous channel estimate h₁ are related according tothe following equation:ŷ=y−F ₁ h ₁ −F ₂ h _(pred).
 5. The method of claim 4 wherein theconjugate gradient algorithm is performed to solve the followingequation:ŷ=Sh wherein h comprises the current channel estimate.
 6. The of claim 4wherein the performing of a conjugate gradient algorithm to determinethe current channel estimate comprises performing the followingalgorithm: (1) ŷ=y−F₁h₁−F₂h_(pred), r₁=Ŝ^(T)ŷ−Ŝ^(T)Ŝh₁ (2) For k=1 to n,iteratively calculate (a) d_(k)=r_(k)+β_(k)d_(k−1) (b)h_(k+1)=h_(k)+α_(k)d_(k) (c) r_(k+1)=r_(k)−α_(k)q_(k−1) where β₁=0,${\beta_{k \geq 2} = \frac{r_{k}^{T} \cdot r_{k}}{r_{k - 1}^{T} \cdot r_{k - 1}}},$where ${\alpha_{k} = \frac{r_{k}^{T} \cdot r_{k}}{d_{k} \cdot q_{k}}},$and where q_(k)=S^(T)Sd_(k).
 7. The method of claim 6 wherein q_(k) isdetermined by forming a first FFT of the matrix Ŝ, by forming a secondFFT of the matrix Ŝ^(T), by forming a third FFT of d_(k), by multiplyingthe first, second, and third FFTs to produce a multiplication result,and by forming an inverse FFT of the multiplication result.
 8. Themethod of claim 6 wherein the forming of a matrix Ŝ from the data scomprises: forming a matrix S from the data s, wherein the matrix Scontains the first, second, and third portions of the data s; and,forming the matrix Ŝ from the matrix S by setting the second and thirdportions of the data s to zero; wherein the forming of a matrix F₁ fromthe data s comprises: forming the matrix F₁ from the matrix S by settingthe first and third portions of the data s to zero; and, wherein theforming of a matrix F₂ from the data s comprises: forming the matrix F₂from the matrix S by setting the first and second portions of the data sto zero.
 9. The method of claim 8 wherein q_(k) is determined by forminga first FFT of the matrix Ŝ, by forming a second FFT of the matrixŜ^(T), by forming a third FFT of d_(k), by multiplying the first,second, and third FFTs to produce a multiplication result, and byforming an inverse FFT of the multiplication result.
 10. The method ofclaim 1 wherein the performing of a conjugate gradient algorithmcomprises determining a quantity q_(k) according to the followingequation:q _(k) =Ŝ ^(T) Ŝd _(k), wherein d_(k) is dependent upon the receivedsignal y, the matrix Ŝ, and the matrices F_(1 and F) ₂, and whereinq_(k) is determined by forming a first FFT of the matrix Ŝ, by forming asecond FFT of the matrix Ŝ^(T), by forming a third FFT of d_(k), bymultiplying the first, second, and third FFTs to produce amultiplication result, and by forming an inverse FFT of themultiplication result.
 11. The method of claim 1 wherein the forming ofa matrix Ŝ from the data s comprises: forming a matrix S from the datas, wherein the matrix S contains the first, second, and third portionsof the data s; and, forming the matrix Ŝ from the matrix S by settingthe second and third portions of the data s to zero; wherein the formingof a matrix F₁ from the data s comprises: forming the matrix F₁ from thematrix S by setting the first and third portions of the data s to zero;and, wherein the forming of a matrix F₂ from the data s comprises:forming the matrix F₂ from the matrix S by setting the first and secondportions of the data s to zero.
 12. The method of claim 11 wherein theperforming of a conjugate gradient algorithm comprises determining aquantity q_(k) according to the following equation:q _(k) =Ŝ ^(T) Ŝd _(k), wherein d_(k) is dependent upon the receivedsignal y, the matrix Ŝ, and the matrices F₁ and F₂, and wherein q_(k) isdetermined by forming a first FFT of the matrix Ŝ, by forming a secondFFT of the matrix Ŝ^(T), by forming a third FFT of d_(k), by multiplyingthe first, second, and third FFTs to produce a multiplication result,and by forming an inverse FFT of the multiplication result.
 13. A methodof processing a received signal y to produce a current channel estimatefrom one or more past channel estimates that include a previous channelestimate h₁ comprising: (a) decoding the received signal y to form datas; (b) forming a convolution matrix Ŝ from a first portion of the datas; (c) forming a matrix F₁ from a second portion the data s, wherein thesecond portion of the data s includes data that is less recent than thedata in the first portion of the data s; (d) forming a matrix F₂ from athird portion the data s, wherein the third portion of the data sincludes data that is more recent than the data in the first portion ofthe data s; (e) determining a predicted channel estimate h_(pred) basedon the one or more past channel estimates; and, (f) performing aconjugate gradient algorithm to determine the current channel estimate,wherein the conjugate gradient algorithm is based on the received signaly, the matrix Ŝ, the matrices F₁ and F₂, the predicted channel estimateh_(pred), and the previous channel estimate h₁, wherein the conjugategradient algorithm includes (i) forming FFTs based on the receivedsignal y, the matrix Ŝ, and the matrices F₁ and F₂, (ii) multiplying theFFTs to form a multiplication product, and (iii) forming an inverse FFTof the multiplication product.
 14. The method of claim 13 wherein thedetermining of a predicted channel estimate h_(pred) comprisesextrapolating the predicted channel estimate from k ones of the pastchannel estimates.
 15. The method of claim 14 wherein k≧2.
 16. Themethod of claim 13 wherein the received signal y, the matrix Ŝ, thematrices F₁ and F₂, the predicted channel estimate h_(pred), and theprevious channel estimate h₁ are related according to the followingequation:ŷ=y−F ₁ h ₁ −F ₂ h _(pred).
 17. The method of claim 16 wherein theconjugate gradient algorithm is performed to solve the followingequation:ŷ=Sh wherein h comprises the current channel estimate.
 18. The method ofclaim 16 wherein the performing of a conjugate gradient algorithm todetermine the current channel estimate comprises performing thefollowing algorithm: (1) ŷ=y−F₁h₁−F₂h_(pred), r₁=Ŝ^(T)ŷ−Ŝ^(T)Ŝh₁ (2) Fork=1 to n, iteratively calculate (a) d_(k)=r_(k)+β_(k)d_(k−1) (b)h_(k+1)=h_(k)+α_(k)d_(k) (c) r_(k+1)=r_(k)−α_(k)q_(k−1) where β₁=0,${\beta_{k \geq 2} = \frac{r_{k}^{T} \cdot r_{k}}{r_{k - 1}^{T} \cdot r_{k - 1}}},$where ${\alpha_{k} = \frac{r_{k}^{T} \cdot r_{k}}{d_{k} \cdot q_{k}}},$where q_(k)=S^(T)Sd_(k).
 19. The method of claim 18 wherein q_(k) isdetermined by forming a first FFT of the matrix Ŝ, by forming a secondFFT of the matrix Ŝ^(T), by forming a third FFT of d_(k), by multiplyingthe first, second, and third FFTs to produce a multiplication result,and by forming an inverse FFT of the multiplication result.
 20. Themethod of claim 18 wherein the forming of a matrix Ŝ from the data scomprises: forming a matrix S from the data s, wherein the matrix Scontains the first, second, and third portions of the data s; and,forming the matrix Ŝ from the matrix S by setting the second and thirdportions of the data s to zero; wherein the forming of a matrix F₁ fromthe data s comprises: forming the matrix F₁ from the matrix S by settingthe first and third portions of the data s to zero; and, wherein theforming of a matrix F₂ from the data s comprises: forming the matrix F₂from the matrix S by setting the first and second portions of the data sto zero.
 21. The method of claim 20 wherein q_(k) is determined byforming a first FFT of the matrix S, by forming a second FFT of thematrix Ŝ^(T), by forming a third FFT of d_(k), by multiplying the first,second, and third FFTs to produce a multiplication result, and byforming an inverse FFT of the multiplication result.
 22. The method ofclaim 13 wherein the performing of a conjugate gradient algorithmcomprises determining a quantity q_(k) according to the followingequation:q _(k) =Ŝ ^(T) Ŝd _(k), wherein d_(k) is dependent upon the receivedsignal y, the matrix Ŝ, and the matrices F_(1 and F) ₂, and whereinq_(k) is determined by forming a first FFT of the matrix Ŝ, by forming asecond FFT of the matrix Ŝ^(T), by forming a third FFT of d_(k), bymultiplying the first, second, and third FFTs to produce amultiplication result, and by forming an inverse FFT of themultiplication result.
 23. The method of claim 13 wherein the forming ofa matrix Ŝ from the data s comprises: forming a matrix S from the datas, wherein the matrix S contains the first, second, and third portionsof the data s; and, forming the matrix Ŝ from the matrix S by settingthe second and third portions of the data s to zero; wherein the formingof a matrix F₁ from the data s comprises: forming the matrix F₁ from thematrix S by setting the first and third portions of the data s to zero;and, wherein the forming of a matrix F₂ from the data s comprises:forming the matrix F₂ from the matrix S by setting the first and secondportions of the data s to zero.
 24. The method of claim 23 wherein theperforming of a conjugate gradient algorithm comprises determining aquantity q_(k) according to the following equation:q _(k) =Ŝ ^(T) Ŝd _(k), wherein d_(k) is dependent upon the receivedsignal y, the matrix Ŝ, and the matrices F₁ and F₂, and wherein q_(k) isdetermined by forming a first FFT of the matrix Ŝ, by forming a secondFFT of the matrix Ŝ^(T), by forming a third FFT of d_(k), by multiplyingthe first, second, and third FFTs to produce a multiplication result,and by forming an inverse FFT of the multiplication result.